Three-dimensional co-prime cubic array direction-of-arrival estimation method based on a cross-correlation tensor

ABSTRACT

The present disclosure discloses a three-dimensional co-prime cubic array direction-of-arrival estimation method based on a cross-correlation tensor, mainly solving the problems of multi-dimensional signal structured information loss and Nyquist mismatch in existing methods and comprising the following implementing steps: constructing a three-dimensional co-prime cubic array; carrying out tensor modeling on a receiving signal of the three-dimensional co-prime cubic array; calculating six-dimensional second-order cross-correlation tensor statistics; deducing a three-dimensional virtual uniform cubic array equivalent signal tensor based on cross-correlation tensor dimension merging transformation; constructing a four-dimensional virtual domain signal tensor based on mirror image augmentation of the three-dimensional virtual uniform cubic array; constructing a signal and noise subspace in a Kronecker product form through virtual domain signal tensor decomposition; and acquiring a direction-of-arrival estimation result based on three-dimensional spatial spectrum search.

TECHNICAL FIELD

The present disclosure belongs to the technical field of array signalprocessing, particularly relates to a statistic signal processingtechnology of a multi-dimensional space sparse array virtual domainstatistics, and more particularly relates to a three-dimensionalco-prime cubic array direction-of-arrival estimation method based on across-correlation tensor, which can be used for target positioning.

BACKGROUND

As a typical systematic sparse array architecture, a co-prime array hasthe advantages of large apertures and high resolution, and can breakthrough the bottleneck of traditional uniform array direction-of-arrivalestimation in an aspect of performance. As array element sparsearrangement of the co-prime array cannot meet the requirement of aNyquist sampling rate, in order to realize Nyquist matcheddirection-of-arrival estimation, it is often to deduct a co-prime arrayreceiving signal to a second-order statistics model to construct anaugmented virtual uniform array and extract angle information by use ofits corresponding virtual domain equivalent signal. To meet requirementsof such fields as radar detection, 5G communication and medical imagingon three-dimensional space target direction finding accuracy, athree-dimensional co-prime cubic array with a larger stereo aperture andits corresponding virtual domain signal processing have attractedextensive attention. In the traditional direction-of-arrival estimationmethod oriented to co-prime linear arrays and co-prime planar arrays, itis often to overlap receiving signals into vectors for processing,deduce a linear virtual domain equivalent signal by vectorizedsecond-order autocorrelation statistics, and introduce spatial smoothingto resolve an identifiable problem of a signal source so as to realizeNyquist matched direction-of-arrival estimation. However, when thisvectorized signal processing method is simply expanded to athree-dimensional co-prime cubic array scene, not only is an originalspatial information structure of a multi-dimensional receiving signaldestroyed, but also the virtual domain signal model deduced by thevectorized signal has the problems of larger linear scale, aliasing andmismatch of multi-dimensional spatial information and the like.

A tensor is a multi-dimensional data type and can be used for storingcomplicated multi-dimensional signal information. Methods includingmulti-dimensional signal characteristic analysis, high-order singularvalue decomposition, tensor decomposition and the like provide richmathematic tools for tensor-oriented signal processing. In recent years,a tensor model has been widely used in array signal processing, imagesignal processing, statistics and other fields. Therefore, by using atensor to construct a three-dimensional co-prime cubic array receivingsignal, the original structure information of the multi-dimensionalreceiving signal can be effectively retained, a Nyquist matched virtualdomain signal processing technology is promoted in a tensor space so asto provide an important theoretical tool for improving performance ofdirection-of-arrival estimation. Whereas the existing tensor method hasnot yet involved in signal processing oriented to three-dimensionalco-prime cubic array virtual domain statistics, and it still adopts atraditional autocorrelation statistics-based calculation method, whichis incapable of effectively matching a traditional linear virtual domaindeduction means into a three-dimensional co-prime cubic array scene.Hence, how to design a multi-dimensional virtual domain model based onstatistics characteristic analysis of the three-dimensional co-primecubic array receiving signal tensor so as to realize Nyquist matcheddirection-of-arrival estimation, becomes a key problem that needs to beurgently resolved at current.

SUMMARY

In view of the problems of multi-dimensional signal structuredinformation loss and Nyquist mismatch in existing methods, an object ofthe present disclosure is to provide a three-dimensional co-prime cubicarray direction-of-arrival estimation method based on across-correlation tensor so as to provide a feasible thought and aneffective solution for establishing relevance between athree-dimensional co-prime cubic array multi-dimensional virtual domainand cross-correlation tensor statistics, and realizing Nyquist matcheddirection-of-arrival estimation by means of structured information ofthe cross-correlation virtual domain signal tensor.

The object of the present disclosure is realized by adopting thefollowing technical solution: a three-dimensional co-prime cubic arraydirection-of-arrival estimation method based on a cross-correlationtensor, comprises the following steps:

(1) constructing at a receiving end with

+

−1 physical antenna array elements in accordance with a structure of athree-dimensional co-prime cubic array, wherein

and

,

and

, and

and

are respectively a pair of co-prime integers, and the three-dimensionalco-prime cubic array is decomposable into two sparse and uniform cubicsubarrays

and

;

(2) supposing there are K far-field narrowband non-coherent signalsources from a direction of {(θ₁, φ₁), (θ₂, φ₂), . . . , (θ_(K),φ_(K))}, carrying out modeling on a receiving signal of the sparse anduniform cubic subarray

of the three-dimensional co-prime cubic array via a four-dimensionaltensor

∈

×

×

×T (T is a number of sampling snapshots) as follows:

=Σ_(k=1) ^(K)

(μ_(k))∘

(ν_(k))∘

(ω_(k))∘s_(k)+

wherein, s_(k)=[s_(k,1), s_(k,2), . . . , s_(k,T)]^(T) is amulti-snapshot sampling signal waveform corresponding to a k^(th)incident signal source, (⋅)^(T) represents a transposition operation, ∘represents an external product of vectors,

∈

×

×

×T is a noise tensor mutually independent from each signal source,

(μ_(k)),

(ν_(k)) and

(ω_(k)) are steering vectors of the three-dimensional sparse and uniformcubic subarray

in the directions of an x axis, a y axis and a z axis respectively, anda signal source corresponding to a direction-of-arrival of (θ_(k),φ_(k))is represented as:

$\begin{matrix}{{{a_{x}^{({\mathbb{Q}}_{i})}( \mu_{k} )} = \lbrack {1,\ e^{{- j}\pi x_{{\mathbb{Q}}_{i}}^{(2)}\mu_{k}},\ldots,\ e^{{- j}\pi x_{{\mathbb{Q}}_{i}}^{(M_{x}^{({\mathbb{Q}}_{i})})}\mu_{k}}} \rbrack^{T}},} & (2)\end{matrix}$a_(y)^((ℚ_(i)))(v_(k)) = [1, e^(−jπy_(ℚ_(i))⁽²⁾v_(k)), …, e^(−jπy_(ℚ_(i))^((M_(y)^((ℚ_(i)))))v_(k))]^(T),a_(z)^((ℚ_(i)))(ω_(k)) = [1, e^(−jπz_(ℚ_(i))⁽²⁾ω_(k)), …, e^(−jπz_(ℚ_(i))^((M_(z)^((ℚ_(i))))))ω_(k)]^(T),

wherein,

(i₁=1,2, . . . ,

),

(i₂=1,2, . . . ,

) and

(i₃=1,2, . . . ,

) respectively represent actual locations of

in i₁ ^(th), i₂ ^(th) and i₃ ^(th) physical antenna array elements inthe directions of the x axis, the y axis and the z axis, and

=

=

=0, μ_(k)=sin(φ_(k))cos(θ_(k)), ν_(k)=sin(φ_(k))sin(θ_(k)),ω_(k)=cos(φ_(k)), j=√{square root over (−1)};

(3) based on four-dimensional receiving signal tensors

and

of the two three-dimensional sparse and uniform cubic subarrays

and

, solving their cross-correlation statistics to obtain a six-dimensionalspace information-covered second-order cross-correlation tensor

×

×

×

×

×

:

${\mathcal{R}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}} = {{E\lbrack {{< \mathcal{X}_{{\mathbb{Q}}_{1}}},{\mathcal{X}_{{\mathbb{Q}}_{2}}^{*} >_{4}}} \rbrack} = {{\sum\limits_{k = 1}^{K}{\sigma_{k}^{2}{{a_{x}^{({\mathbb{Q}}_{1})}( \mu_{k} )} \circ {a_{y}^{({\mathbb{Q}}_{1})}( v_{k} )} \circ {a_{z}^{({\mathbb{Q}}_{1})}( \omega_{k} )} \circ {a_{x}^{{({\mathbb{Q}}_{2})}^{*}}( \mu_{k} )} \circ {a_{y}^{{({\mathbb{Q}}_{2})}^{*}}( v_{k} )} \circ {a_{z}^{{({\mathbb{Q}}_{2})}^{*}}( \omega_{k} )}}}} + \mathcal{N}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}}}}},$

wherein, σ_(k) ²=E[s_(k)s_(k)*] represents a power of a k^(th) incidentsignal source,

=E[<

,

>₄] represents a six-dimensional cross-correlation noise tensor,<⋅,⋅>_(r) represents a tensor contraction operation of two tensors alonga r^(th) dimension, E[⋅] represents an operation of taking a mathematicexpectation, and (⋅)* represents a conjugate operation. Asix-dimensional tensor

merely has an element with a value of σ_(n) ² in a (1,1,1,1,1,1)^(th)location, σ_(n) ² representing a noise power, and with a value of 0 inother locations;

(4) as a first dimension and a fourth dimension (represented by steeringvectors

(μ_(k)) and

(μ_(k))) of the cross-correlation tensor

represent space information in the direction of the x axis, a seconddimension and a fifth dimension (represented by steering vectors

(ν_(k)) and

(ν_(k))) represent space information in the direction of the y axis, anda third dimension and a sixth dimension (represented by steering vectors

(ω_(k)) and

(ω_(k))) represent space information in the direction of the z axis,defining dimension sets

={1,4}

={2,5} and

={3,6}, and carrying out tensor transformation of dimension merging onthe cross-correlation tensor

to obtain a virtual domain second-order equivalent signal tensor

∈

×

×

:

${\mathcal{U}_{\mathbb{W}}\overset{\Delta}{=}{\mathcal{R}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{{2\lbrack{{\mathbb{J}}_{1},{\mathbb{J}}_{2},{\mathbb{J}}_{3}}}\}}} = {\sum_{k = 1}^{K}{\sigma_{k}^{2}{{b_{x}( \mu_{k} )} \circ {b_{y}( v_{k} )} \circ {b_{z}( \omega_{k} )}}}}}},$

wherein, b_(x)(μ_(k))=

(μ_(k))⊗

(μ_(k)), b_(y)(ν_(k))=

(ν_(k))⊗

(ν_(k)) and b_(z)(ω_(k))=

(ω_(k))⊗

(ω_(k)) respectively construct augmented virtual arrays in thedirections of the x axis, the y axis and the z axis through formingarrays of difference sets on exponential terms, b_(x)(μ_(k)),b_(y)(ν_(k)) and b_(z)(ω_(k)) are respectively equivalent to steeringvectors of the virtual arrays in the x axis, the y axis and the z axisto correspond to signal sources in a direction-of-arrival of(θ_(k),φ_(k)), and ⊗ represents a product of Kronecker. Therefore,

corresponds to an augmented three-dimensional virtual non-uniform cubicarray

; to simplify a deduction process, the six-dimensional noise tensor

is omitted in a theoretical modeling step about

;

comprises a three-dimensional uniform cubic array

with (3

−

+1)×(3

−

+1)×(3

−

+1) virtual array elements, represented as:

={(x,y,z)|x=p _(x) d,y=p _(y) d,z=p _(z) d,−

≤p _(x)≤−

+2

, −

≤p _(y)≤−

+2

,−

≤p _(z)≤−

+2

};

the equivalent signal tensor

∈

−

+1)×(3

−

+1)×(3

−

+1) of the three-dimensional uniform cubic array

is modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘ b _(y)(ν_(k))∘ b _(z)(ω_(k)),

wherein,

${{{\overset{\_}{b}}_{x}( \mu_{k} )} = {\lbrack {e^{{- j}{\pi({- M_{x}^{({\mathbb{Q}}_{1})}})}\mu_{k}},e^{{- j}{\pi({{- M_{x}^{({\mathbb{Q}}_{1})}} + 1})}\mu_{k}},\ldots,e^{{- j}{\pi({{- M_{x}^{({\mathbb{Q}}_{2})}} + {2M_{x}^{({\mathbb{Q}}_{1})}}})}\mu_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{x}^{({\mathbb{Q}}_{1})}} - M_{x}^{({\mathbb{Q}}_{2})} + 1}}},$${{\overset{\_}{b}}_{y}( v_{k} )} = {\lbrack {e^{{- j}{\pi({- M_{y}^{({\mathbb{Q}}_{1})}})}v_{k}},e^{{- j}{\pi({{- M_{y}^{({\mathbb{Q}}_{1})}} + 1})}v_{k}},\ldots,e^{{- j}{\pi({{- M_{y}^{({\mathbb{Q}}_{2})}} + {2M_{y}^{({\mathbb{Q}}_{1})}}})}v_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{y}^{({\mathbb{Q}}_{1})}} - M_{y}^{({\mathbb{Q}}_{2})} + 1}}$and${{\overset{\_}{b}}_{z}( \omega_{k} )} = {\lbrack {e^{{- j}{\pi({- M_{z}^{({\mathbb{Q}}_{1})}})}\omega_{k}},e^{{- j}{\pi({{- M_{z}^{({\mathbb{Q}}_{1})}} + 1})}\omega_{k}},\ldots,e^{{- j}{\pi({{- M_{z}^{({\mathbb{Q}}_{2})}} + {2M_{z}^{({\mathbb{Q}}_{1})}}})}\omega_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{z}^{({\mathbb{Q}}_{1})}} - M_{z}^{({\mathbb{Q}}_{2})} + 1}}$

respectively represent steering vectors of the three-dimensional virtualuniform cubic array

in the x axis, the y axis and the z axis corresponding to signal sourcesin the direction-of-arrival of (θ_(k),φ_(k));

(5) as a mirror image portion

_(sym) of the three-dimensional virtual uniform cubic array

is represented as:

_(sym)={(x,y,z)|x={hacek over (p)} _(x) d,y={hacek over (p)} _(y)d,z={hacek over (p)} _(z) d,

−2

≤{hacek over (p)} _(x)≤

,

−2

≤{hacek over (p)} _(y)≤

,

−2

≤{hacek over (p)} _(z)≤

};

carrying out transformation by using the equivalent signal tensor

of the three-dimensional virtual uniform cubic array W to obtain anequivalent signal tensor

_(sym)∈

−

+1)×(3

−

+1)×(3

−

+1) of a three-dimensional mirror image virtual uniform cubic array

_(sym), specifically comprising: carrying out a conjugate operation onthe three-dimensional virtual domain signal tensor

to obtain

, carrying out position reversal on elements in the

along directions of three dimensions successively so as to obtain theequivalent signal tensor

_(sym) corresponding to the

_(sym);

superposing the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

and the equivalent signal tensor

_(sym) of the mirror image virtual uniform cubic array

_(sym) in the fourth dimension (i.e., a dimension representing mirrorimage transformation information) to obtain a four-dimensional virtualdomain signal tensor

∈

−

+1)×(3

−

+1)×(1

−

+1)×2, modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘ b _(y)(ν_(k))∘ b_(z)(ω_(k))∘c(μ_(k), ν_(k), ω_(k)),

wherein,

c(μ_(k), v_(k), ω_(k)) = [1, e^(−jπ((M_(x)^((ℚ₂)) − M_(x)^((ℚ₁)))μ_(k) + (M_(y)^((ℚ₂)) − M_(y)^((ℚ₁)))v_(k) + (M_(z)^((ℚ₂)) − M_(z)^((ℚ₁)))ω_(k)))]^(T)

is a three-dimensional space mirror image transformation factor vector;

(6) carrying out CANDECOMP/PARACFAC decomposition on thefour-dimensional virtual domain signal tensor

to obtain factor vectors b _(x)(μ_(k)), b _(y)(ν_(k)), b _(z)(ω_(k)) andc(μ_(k),ν_(k),ω_(k)), k=1,2, . . . , K, corresponding tofour-dimensional space information, and constructing a signal subspaceV_(s)∈

^(V×K) through a form of their Kronecker products:

V _(s)=orth([ b _(x)(μ₁)⊗ b _(y)(ν₁)⊗ b _(z)(ω₁)⊗c(μ₁, ν₁, ω₁), b_(x)(μ₂)⊗ b _(y)(ν₂)⊗ b _(z)(ω₂)⊗c(μ₂, ν₂, ω₂), . . . , b _(x)(μ_(K))⊗ b_(y)(ν_(K))⊗ b _(z)(ω_(K))⊗c(μ_(K), ν_(K), ω_(K))]),

wherein, orth(⋅) represents a matrix orthogonalization operation, V=2(3

−

+1)(3

−

+1)(3

−

+1); by using V_(n)∈

^(V×(V−K)) to represent a noise subspace, V_(n)V_(n) ^(H) is obtained byV_(s):

V _(n) V _(n) ^(H) =I−V _(s) V _(s) ^(H),

wherein, I represents a unit matrix; (⋅)^(H) represents a conjugatetransposition operation; and

(7) traversing a two-dimensional direction-of-arrival of ({tilde over(θ)},{tilde over (φ)}), calculating corresponding parameters {tilde over(ν)}_(k)=sin({tilde over (φ)}_(k))cos({tilde over (θ)}_(k)), {tilde over(ν)}_(k)=sin({tilde over (φ)}_(k))sin({tilde over (θ)}_(k)) and {tildeover (ω)}_(k)=cos ({tilde over (φ)}_(k)), and constructing a steeringvector

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))∈

^(V) corresponding to the three-dimensional virtual uniform cubic array

, represented as:

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))= b_(x)({tilde over (μ)}_(k))⊗ b _(y)({tilde over (ν)}_(k))⊗ b _(z)({tildeover (ω)}_(k))⊗c({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over(ω)}_(k)),

wherein, {tilde over (θ)}∈[−90°, 90°], {tilde over (φ)}∈[0°, 180°]. Athree-dimensional spatial spectrum

({tilde over (θ)},{tilde over (ω)}) is calculated as follows:

({tilde over (θ)},{tilde over (φ)})=1/(

^(H)({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))(V_(n) V _(n) ^(H))

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))).

Spectral peak search is carried out on the three-dimensional spatialspectrum

({tilde over (θ)},{tilde over (φ)}) to obtain a direction-of-arrivalestimation result.

Further, the structure of the three-dimensional co-prime cubic array instep (1) is specifically described as: a pair of three-dimensionalsparse and uniform cubic subarrays

and

are constructed in a rectangular coordinate system, wherein

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, the y axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(1x),

dm_(1y),

dm_(1z)), m_(1x)=0,1, . . . ,

−1, m_(1y)=0, 1, . . . ,

−1, m_(1z)=0,1, . . . ,

−1};

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, the y axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(2x),

dm_(2y),

dm_(2z)), m_(2x)=0,1, . . . ,

−1, m_(2y)=0,1, . . . ,

−1, m_(2z)=0,1, . . . ,

−1}; a unit spacing d has a value half of an incident narrowband signalwavelength λ, i.e., d=λ/2 subarray combination is carried out on the

and

in such a way that array elements on the (0,0,0) location in therectangular coordinate system are overlapped so as to obtain athree-dimensional co-prime cubic array actually containing

+

−1 physical antenna array elements.

Further, the second-order cross-correlation tensor statistics of thethree-dimensional co-prime cubic array in step (3) are estimated bycalculating cross-correlation statistics of T sampling snapshots of thereceiving signal tensors

(t) and

(t) in reality:

${\overset{\hat{}}{\mathcal{R}}}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}} = {\frac{1}{T}{\sum_{t = 1}^{T}{{{\mathcal{X}_{{\mathbb{Q}}_{1}}(t)} \circ {\mathcal{X}_{{\mathbb{Q}}_{2}}^{*}(t)}}.}}}$

Further, the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

in step (4) can be obtained by selecting elements in the equivalentsignal tensor

of the three-dimensional virtual non-uniform cubic array

corresponding to locations of virtual array elements in the

.

Further, in step (6), CANDECOMP/PARAFAC decomposition is carried out onthe four-dimensional virtual domain signal tensor

to obtain factor matrixes B _(x)=[b _(x)(μ₁), b _(x)(μ₂), . . . , b_(x)(μ_(K))], B _(y)=[b _(y)(ν₁), b _(y)(ν₂), . . . , b _(y)(ν_(K))], B_(z)=[b _(z)(ω₁), b _(z)(ω₂), . . . , b _(z)(ω_(K))] and C=[c(μ₁,ν₁,ω₁),c(μ₂,ν₂,ω₂), . . . , c(μ_(K),ν_(K),ω_(K))], wherein, CANDECOMP/PARAFACdecomposition of the four-dimensional virtual domain signal tensor

follows a uniqueness condition as follows:

_(rank)( B _(x))+

_(rank)( B _(y))+

_(rank)( B _(z))+

_(rank)(C)≥2K+3,

wherein,

_(rank)(⋅) represents a Kruskal rank of a matrix, and

_(rank)(B _(x))=min (3

−

+1, K),

_(rank)(B _(y))=min (3

−

+1, K),

_(rank)(B _(z))=min (3

−

−1, K),

_(rank)(C)=min (2, K), min (⋅) represents an operation of taking aminimum value; when spatial smoothing is not introduced to process thededuced four-dimensional virtual domain signal tensor

, a uniqueness inequation of the above CANDECOMP/PARACFAC decompositionis established, indicating that angle information of a signal source canbe effectively extracted in the method of the present disclosure in noneed of a spatial smoothing step.

Further, in step (7), a process of obtaining a direction-of-arrivalestimation result by three-dimensional spatial spectrum searchspecifically comprises: fixing a value of {tilde over (φ)} at 0°,gradually increasing {tilde over (θ)} to 90° from −90° at an interval of0.1°, increasing the {tilde over (φ)} to 0.1° from 0°, increasing the{tilde over (θ)} to 90° from −90° at an interval of 0.1° once again, andrepeating this process until the {tilde over (φ)} is increased to 180°,calculating a corresponding

({tilde over (θ)},{tilde over (φ)}) in each two-dimensionaldirection-of-arrival of ({tilde over (θ)},{tilde over (φ)}) so as toconstruct a three-dimensional spatial spectrum in a two-dimensionaldirection-of-arrival plane; and searching peak values of thethree-dimensional spatial spectrum

({tilde over (θ)},{tilde over (φ)}) in the two-dimensionaldirection-of-arrival plane, permutating response values corresponding tothese peak values in a descending order, and taking two-dimensionalangle values corresponding to first K spectral peaks as thedirection-of-arrival estimation result of a corresponding signal source.

Compared with the prior art, the present disclosure has the followingadvantages:

(1) a tensor is used to represent a multi-dimensional receiving signalof a three-dimensional co-prime cubic array in the present disclosure,which, compared with a traditional vectorized signal processing method,effectively retains the original structured information ofmulti-dimensional receiving signals, de-structures space information andtime information of a tensor signal using a tensor algebra tool, andavoids aliasing of multi-dimensional time-space information;

(2) in the present disclosure, a multi-dimensional virtual domainequivalent signal tensor is deduced based on six-dimensionalcross-correlation tensor statistics, solving the failure problem of atraditional autocorrelation statistics-based linear virtual domaindeduction method in a three-dimensional co-prime cubic array scene,establishing relevance between co-prime tensor signals andmulti-dimensional virtual domains, and laying a foundation for realizingNyquist matched direction-of-arrival estimation;

(3) in the present disclosure, usable information of thethree-dimensional co-prime cubic array virtual domain is fully mined, adirectly decomposable four-dimensional virtual domain signal tensor isconstructed by mirror image augmentation of the three-dimensionalvirtual uniform cubic array, and direction-of-arrival estimation of asignal source is realized without introducing a spatial smoothing step.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is an overall flow block diagram of the present disclosure.

FIG. 2 is a schematically structural diagram of a three-dimensionalco-prime cubic array of the present disclosure.

FIG. 3 is a schematically structural diagram of a deducedthree-dimensional virtual uniform cubic array of the present disclosure.

FIG. 4 is an effect diagram of direction-of-arrival estimation of atraditional method based on vectorized signal processing.

FIG. 5 is an effect diagram of direction-of-arrival estimation of amethod put forward by the present disclosure.

DETAILED DESCRIPTION OF THE EMBODIMENTS

The technical solution of the present disclosure will be furtherexplained in detail by referring to the appended drawings.

In order to solve the problems of multi-dimensional signal structuredinformation loss and Nyquist mismatch in existing methods, the presentdisclosure puts forwards a three-dimensional co-prime cubic arraydirection-of-arrival estimation method based on a cross-correlationtensor. In combination with cross-correlation tensor statistic analysis,multi-dimensional virtual domain tensor space extension,cross-correlation virtual domain signal tensor decomposition and othermeans, relevance between three-dimensional co-prime cubic arraycross-correlation tensor statistics and a virtual domain is establishedto realize Nyquist matched two-dimensional direction-of-arrivalestimation. As shown in FIG. 1 , the present disclosure comprises thefollowing implementing steps:

step 1: constructing a three-dimensional co-prime cubic array. Thethree-dimensional co-prime cubic array is constructed at a receiving endwith

+

−1 physical antenna array elements; as shown in FIG. 2 , a pair ofthree-dimensional sparse and uniform cubic subarrays

and

are constructed in a rectangular coordinate system, wherein

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, the y axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(1x),

dm_(1y),

dm_(1z)), m_(1x)=0,1, . . . ,

−1, m_(1y)=0,1, . . . ,

−1, m_(1z)=0,1, . . . ,

−1};

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, the y axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(2x),

dm_(2y),

dm_(2z)), m_(2x)=0,1, . . . ,

−1, m_(2y)=0,1, . . . ,

−1, 3_(2z)=0,1, . . . ,

−1}; wherein

and

,

and

, and

and

are respectively a pair of co-prime integers; a unit spacing d has avalue half of an incident narrowband signal wavelength λ, i.e., d=λ/2;subarray combination is carried out on the

and

in such a way that array elements on the (0,0,0) location in therectangular coordinate system are overlapped so as to obtain athree-dimensional co-prime cubic array actually containing

+

−1 physical antenna array elements.

step 2: carrying out tensor modeling on a receiving signal of thethree-dimensional co-prime cubic array. Supposing there are K far-fieldnarrowband non-coherent signal sources from a direction of {(θ₁,φ₁),(θ₂,φ₂), . . . , (θ_(K),φ_(K))}, a sampling snapshot signal at the timet of the sparse and uniform cubic subarray

of the three-dimensional co-prime cubic array is represented by athree-dimensional space information-covered tensor

(t)∈

×

×

, receiving signal tensors

(t) of T sampling snapshots are superposed in a fourth dimension (i.e.,time dimension) to obtain a four-dimensional receiving signal tensor

∈

×

×

×T corresponding to the sparse and uniform cubic subarray

, modeled as:

=Σ_(k=1) ^(K)

(μ_(k))∘

(ν_(k))∘

(ω_(k))∘s _(k)+

,

wherein, s_(k)=[s_(k,1), s_(k,2), . . . , s_(k,T)]^(T) is amulti-snapshot sampling signal waveform corresponding to a k^(th)incident signal source, (⋅)^(T) represents a transposition operation, ∘represents an external product of vectors,

×

×

×T is a noise tensor mutually independent from each signal source,

(μ_(k)),

(ν_(k)) and

(ω_(k)) are steering vectors of the three-dimensional sparse and uniformcubic subarray

in the directions of the x axis, the y axis and the z axis respectively,and a signal source corresponding to a direction-of-arrival of(θ_(k),φ_(k)) is represented as:

a_(x)^((ℚ_(i)))(μ_(k)) = [1, e^(−jπx_(ℚ_(i))⁽²⁾μ_(k)), …, e^(−jπx_(ℚ_(i))^((M_(x)^((ℚ_(i)))))μ_(k))]^(T),a_(y)^((ℚ_(i)))(v_(k)) = [1, e^(−jπy_(ℚ_(i))⁽²⁾v_(k)), …, e^(−jπy_(ℚ_(i))^((M_(y)^((ℚ_(i)))))v_(k))]^(T),a_(z)^((ℚ_(i)))(ω_(k)) = [1, e^(−jπz_(ℚ_(i))⁽²⁾ω_(k)), …, e^(−jπz_(ℚ_(i))^((M_(z)^((ℚ_(i))))))ω_(k)]^(T),

wherein,

(i₁=1,2, . . . ,

) ,

(i₂=1,2, . . . ,

) and

(i₃=1,2, . . . ,

) respectively represent actual locations of

in i₁ ^(th), i₂ ^(th) and i₃ ^(th) physical antenna array elements inthe directions of the x axis, the y axis and the z axis, and

=

=

=0, μ_(k)=sin(φ_(k))cos(θ_(k)), ν_(k)=sin(φ_(k))sin(θ_(k)),ω_(k)=cos(φ_(k)), j=√{square root over (−1)};

step 3: calculating six-dimensional second-order cross-correlationtensor statistics. As receiving signal tensors

and

of the two subarrays

and

meet characteristics of co-prime numbers in structure size, it is unableto superpose

and

into a tensor signal and then to calculate its second-orderautocorrelation statistics. Therefore, their cross-correlationstatistics are solved to obtain a six-dimensional spaceinformation-covered second-order cross-correlation tensor

∈

×

×

×

×

×

, calculated as:

ℛ_(ℚ₁ℚ₂) = E[ < 𝒳_(ℚ₁), 𝒳_(ℚ₂)^(*)>₄]

wherein, σ_(k) ²=E[s_(k)s_(k)*] represents a power of a k^(th) incidentsignal source,

=E[<

,

>₄] represents a six-dimensional cross-correlation noise tensor,<⋅,⋅>_(r) represents a tensor contraction operation of two tensors alonga r^(th) dimension, E[⋅] represents an operation of taking a mathematicexpectation, and (⋅)* represents a conjugate operation. Asix-dimensional tensor

merely has an element with a value of σ_(n) ² on a (1,1,1,1,1,1)^(th)location, σ_(n) ² representing a noise power, and with a value of 0 onother locations. In fact, their cross-correlation statistics ofreceiving signal tensors

(t) and

(t) of T sampling snapshots are solved to obtain a second-order samplingcross-correlation tensor

∈

×

×

×

×

×

:

${{\hat{\mathcal{R}}}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}} = {\frac{1}{T}{\sum_{t = 1}^{T}{{\mathcal{X}_{{\mathbb{Q}}_{1}}(t)} \circ {\mathcal{X}_{{\mathbb{Q}}_{2}}^{*}(t)}}}}};$

step 4: deducing a three-dimensional virtual uniform cubic arrayequivalent signal tensor based on cross-correlation tensor dimensionmerging transformation. As a cross-correlation tensor

contains three-dimensional space information respectively correspondingto the two sparse and uniform cubic subarrays

and

, and by merging dimensions representing space information of a samedirection in the

, arrays of difference sets are formed on exponential terms of thesteering vectors corresponding to the two co-prime subarrays so as toconstruct an augmented virtual array in a three-dimensional space.Specifically, as a first dimension and a fourth dimension (representedby steering vectors

(μ_(k)) and

(μ_(k))) of the cross-correlation tensor

represent space information in the direction of the x axis, a seconddimension and a fifth dimension (represented by steering vectors

(ν_(k)) and

(ν_(k))) represent space information in the direction of the y axis, anda third dimension and a sixth dimension (represented by steering vectors

(ω_(k)) and

(ω_(k))) represent space information in the direction of the z axis,dimension sets

={1,4},

={2,5} and

={3,6} are defined, and tensor transformation of dimension merging iscarried out on the cross-correlation tensor

to obtain a virtual domain second-order equivalent signal tensor

∈

×

×

:

${\mathcal{U}_{\mathbb{W}}\overset{\Delta}{=}{\mathcal{R}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{{2\lbrack{{\mathbb{J}}_{1},{\mathbb{J}}_{2},{\mathbb{J}}_{3}}}\}}} = {\sum_{k = 1}^{K}{\sigma_{k}^{2}{{b_{x}( \mu_{k} )} \circ {b_{y}( v_{k} )} \circ {b_{z}( \omega_{k} )}}}}}},$

wherein, b_(x)(μ_(k))=

(μ_(k))⊗

(μ_(k)), b_(y)(ν_(k))=

(ν_(k))⊗

(ν_(k)) and b_(z)(ω_(k))=

(ω_(k))⊗

(ω_(k)) respectively construct augmented virtual arrays in thedirections of the x axis, the y axis and the z axis through formingarrays of difference sets on exponential terms, b_(x)(μ_(k)),b_(y)(ν_(k)) and b_(z)(ω_(k)) are respectively equivalent to steeringvectors of the virtual arrays in the x axis, the y axis and the z axisto correspond to signal sources in a direction-of-arrival of(θ_(k),φ_(k)) and ⊗ represents a product of Kronecker. Therefore,

corresponds to an augmented three-dimensional virtual non-uniform cubicarray

. To simplify a deduction process, a six-dimensional noise tensor

is omitted in a theoretical modeling step about

. However, in fact, by replacing theoretical cross-correlation tensorstatistics

with sampling cross-correlation tensor statistics

,

is naturally covered in a cross-correlation tensor signal statisticprocessing process;

comprises a three-dimensional uniform cubic array

with (3

−

+1)×(3

−

+1)×(3

−

+1) virtual array elements, represented as:

={(x,y,z)|x=p _(x) d,y=p _(y) d,z=p _(z) d,−

≤p _(x)≤−

+2

, −

≤p _(y)≤−

+2

, −

≤p _(z)≤−

+2

}.

The equivalent signal tensor

∈

−

+1)×(3

−

+1)×(3

−

+1) of the three-dimensional virtual uniform cubic array

can be obtained by selecting elements in

that correspond to locations of virtual array elements in the

, modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘ b _(y)(ν_(k))∘ b _(z)(ω_(k)),

wherein,

${{{\overset{\_}{b}}_{x}( \mu_{k} )} = {\lbrack {e^{- j{\pi({- M_{x}^{({\mathbb{Q}}_{1})}})}\mu_{k}},e^{- j{\pi({{- M_{x}^{({\mathbb{Q}}_{1})}} + 1})}\mu_{k}},\ldots,e^{- j{\pi({{- M_{x}^{({\mathbb{Q}}_{2})}} + {2M_{x}^{({\mathbb{Q}}_{1})}}})}\mu_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{x}^{({\mathbb{Q}}_{1})} - M_{x}^{({\mathbb{Q}}_{2})}} + 1}}},{{{\overset{\_}{b}}_{y}( v_{k} )} = {{\lbrack {e^{- j{\pi({- M_{y}^{({\mathbb{Q}}_{1})}})}v_{k}},e^{- j{\pi({{- M_{y}^{({\mathbb{Q}}_{1})}} + 1})}v_{k}},\ldots,e^{- j{\pi({{- M_{y}^{({\mathbb{Q}}_{2})}} + {2M_{y}^{({\mathbb{Q}}_{1})}}})}v_{k}}} \rbrack^{T} \in {{\mathbb{C}}^{{3M_{y}^{({\mathbb{Q}}_{1})}} - M_{y}^{({\mathbb{Q}}_{2})} + 1}{and}{{\overset{\_}{b}}_{z}( \omega_{k} )}}} = {\lbrack {e^{- j{\pi({- M_{z}^{({\mathbb{Q}}_{1})}})}\omega_{k}},e^{- j{\pi({{- M_{z}^{({\mathbb{Q}}_{1})}} + 1})}\omega_{k}},\ldots,e^{- j{\pi({{- M_{z}^{({\mathbb{Q}}_{2})}} + {2M_{z}^{({\mathbb{Q}}_{1})}}})}\omega_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{z}^{({\mathbb{Q}}_{1})}} - M_{z}^{({\mathbb{Q}}_{2})} + 1}}}}$

respectively represent steering vectors of the three-dimensional virtualuniform cubic array

in the x axis, the y axis and the z axis corresponding to signal sourcesin the direction-of-arrival of (θ_(k),φ_(k));

step 5: constructing a four-dimensional virtual domain signal tensorbased on mirror image augmentation of the three-dimensional virtualuniform cubic array. As the three-dimensional virtual uniform cubicarray

obtained based on cross-correlation tensor dimension mergingtransformation is not symmetric about a coordinate axis, in order toincrease an effective aperture of a virtual array, considering that amirror image portion

_(sym) of the three-dimensional virtual uniform cubic array

is represented as:

_(sym)={(x,y,z)|x={hacek over (p)} _(x) d,y={hacek over (p)} _(y)d,z={hacek over (p)} _(z) d,

−2

≤{hacek over (p)} _(x)≤

,

−2

≤{hacek over (p)} _(y)≤

,

−2

≤{hacek over (p)} _(z)≤

}.

Transformation is carried out by using the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

to obtain an equivalent signal tensor

_(sym)∈

−

+1)×(3

−

+1)×(3

−

+1) of a mirror image virtual uniform cubic array

_(sym). It specifically comprises: a conjugate operation is carried outon the three-dimensional virtual domain signal tensor

to obtain

, position reversal is carried out on elements in the

along directions of three dimensions successively so as to obtain theequivalent signal tensor

_(sym) corresponding to the

_(sym);

the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

and the equivalent signal tensor

_(sym) of the mirror image virtual uniform cubic array

_(sym) are superposed in the fourth dimension (i.e., a dimensionrepresenting mirror image transformation information) to obtain afour-dimensional virtual domain signal tensor

∈

−

+1)×(3

−

+1)×(3

−

+1)×2, modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘ b _(y)(ν_(k))∘ b_(z)(ω_(k))∘c(μ_(k),ν_(k),ω_(k)),

wherein,

c(μ_(k), v_(k), ω_(k)) = [1, e^(−jπ((M_(x)^((ℚ₂)) − M_(x)^((ℚ₁)))μ_(k) + (M_(y)^((ℚ₂)) − M_(y)^((ℚ₁)))v_(k) + (M_(z)^((ℚ₂)) − M_(z)^((ℚ₁)))ω_(k)))]^(T)

is a three-dimensional space mirror image transformation factor vector;

step 6: constructing a signal and noise subspace in the form of aKronecker product by virtual domain signal tensor decomposition.CANDECOMP/PARACFAC decomposition is carried out on the four-dimensionalvirtual domain signal tensor

to obtain factor vectors b _(x)(μ_(k)), b _(y)(ν_(k)), b _(z)(ω_(k)) andc(μ_(k),ν_(k),ω_(k)), k=1,2, . . . , K, corresponding tofour-dimensional space information, and B _(x)=[b _(x)(μ₁) b _(x)(μ₂), .. . , b _(x)(μ_(K))], B _(y)=[b _(y)(ν₁), b _(y)(ν₂), . . . , b_(y)(ν_(K))],B _(z)=[b _(z)(ω₁), b _(z)(ω₂), . . . , b _(z)(ω_(K))] andC=[c(μ₁,ν₁,ω₁), c(μ₂,ν₂ω₂), . . . , c(μ_(K), ν_(K),ω_(K))] are used torepresent factor matrixes. At this point, CANDECOMP/PARAFACdecomposition of the four-dimensional virtual domain signal tensor

follows a uniqueness condition as follows:

_(rank)( B _(x))+

_(rank)( B _(y))+

_(rank)( B _(z))+

_(rank)(C)≥2K+3,

wherein,

_(rank)(⋅) represents a Kruskal rank of a matrix, and

_(rank)(B _(x))=min (3

−

+1, K),

_(rank)(B _(y))=min (3

−

+1, K),

_(rank)(B _(z))=min(3

−

+1, K),

_(rank)(C)=min(2, K), min (⋅) represents an operation of taking aminimum value; when spatial smoothing is not introduced to process thededuced four-dimensional virtual domain signal tensor

, a uniqueness inequation of the above CANDECOMP/PARACFAC decompositionis established, indicating that angle information of a signal source canbe effectively extracted by the method in the present disclosure in noneed of a spatial smoothing step. Further, factor vectors b _(x)(μ_(k)),b _(y)(ν_(k)),b _(z)(ω_(k)) and c(μ_(k),ν_(k),ω_(k)) are obtained bytensor decomposition, and a signal subspace V_(s)∈

^(V×K) is constructed through a form of their Kronecker products:

V _(s)=orth([ b _(x)(μ₁)⊗ b _(y)(ν₁)⊗b _(z)(ω₁)⊗c(μ₁,ν₁,ω₁), b _(x)(μ₂)⊗b _(y)(ν₂)⊗ b _(z)(ω₂)⊗c(μ₂,ν₂,ω₂), . . . , b _(x)(μ_(K))⊗ b_(y)(ν_(K))⊗ b _(z)(ω_(K))⊗c(μ_(K),ν_(K),ω_(K))]),

wherein, orth(⋅) represents a matrix orthogonalization operation, V=2(3

−

+1)(3

−

+1)(3

−

+1); by using V_(n)∈

^(V×(V−K)) to represent a noise subspace, V_(n)V_(n) ^(H) is obtained byV_(s):

V _(n) V _(n) ^(H) =I−V _(s) V _(s) ^(H),

wherein, I represents a unit matrix; (⋅)^(H) represents a conjugatetransposition operation;

step 7: obtaining a direction-of-arrival estimation result based onthree-dimensional spatial spectrum search. A two-dimensionaldirection-of-arrival of ({tilde over (θ)},{tilde over (φ)}) istraversed, corresponding parameters {tilde over (μ)}_(k)=sin({tilde over(φ)}_(k))cos({tilde over (θ)}_(k)), {tilde over (ν)}_(k)=sin({tilde over(φ)}_(k))sin({tilde over (θ)}_(k)) and {tilde over (ω)}_(k)=cos({tildeover (φ)}_(k)) are calculated, and a steering vector

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))∈

^(V) corresponding to the three-dimensional virtual uniform cubic arrayW is constructed, represented as:

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))= b_(x)(μ_(k))⊗ b _(y)({tilde over (ν)}_(k))⊗{tilde over (b)} _(z)({tildeover (ω)}_(k))⊗c({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over(ω)}_(k)),

wherein, {tilde over (θ)}∈[−90°,90°], {tilde over (φ)}∈[0°,180°]. Athree-dimensional spatial spectrum

({tilde over (θ)},{tilde over (φ)}) is calculated as follows:

({tilde over (θ)},{tilde over (φ)})=1/(

^(H)({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))(V_(n) V _(n) ^(H))

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))).

Spectral peak search is carried out on the three-dimensional spatialspectrum

({tilde over (θ)},{tilde over (φ)}) to obtain a direction-of-arrivalestimation result. It specially comprises: a value of {tilde over (φ)}is fixed at 0°, {tilde over (θ)} is gradually increased to 90° from −90°at an interval of 0.1°, then the {tilde over (φ)} is increased to 0.1°from 0°, the {tilde over (θ)} is increased to 90° from −90° at aninterval of 0.1° once again, and this process is repeated until the{tilde over (φ)} is increased to 180°, a corresponding

({tilde over (θ)},{tilde over (φ)}) is calculated in eachtwo-dimensional direction-of-arrival of ({tilde over (θ)},{tilde over(φ)}) so as to construct a three-dimensional spatial spectrum in atwo-dimensional direction-of-arrival plane; peak values of thethree-dimensional spatial spectrum

({tilde over (θ)},{tilde over (φ)}) are searched in the two-dimensionaldirection-of-arrival plane, response values corresponding to these peakvalues are permutated in a descending order, and two-dimensional anglevalues corresponding to first K spectral peaks are taken as thedirection-of-arrival estimation result of a corresponding signal source.

The effects of the present disclosure will be further described withconjunction of a simulation example.

In the simulation example, a three-dimensional co-prime cubic array isused to receive an incident signal, with selected parameters being

=

=

=2,

=

=

=3, i.e., the constructed three-dimensional co-prime cubic arraycomprises

+

−1=34 physical array elements in total. Supposing there are two incidentnarrowband signals, azimuths and pitch angles in the incident directionare respectively [25°,20°] and [45°,40°], and when SNR=0 dB, 800sampling snapshots are adopted for a simulation experiment.

An estimation result of a traditional three-dimensional co-prime cubicarray direction-of-arrival estimation method based on vectorized signalprocessing is as shown in FIG. 4 , wherein the x axis and the y axisrepresent a pitch angle and an azimuth of an incident signal sourcerespectively. An estimation result of a three-dimensional co-prime cubicarray direction-of-arrival estimation method based on across-correlation tensor provided by the present disclosure is as shownin FIG. 5 . It is seen upon comparison that, the method of the presentdisclosure can preciously estimate the two incident signal sources,while the traditional method based on vectorized signal processing isunable to effectively identify the two incident signal sources, showingthe advantages of the method of the present disclosure indirection-of-arrival estimation on resolution and performance.

In conclusion, in the present disclosure, relevance between athree-dimensional co-prime cubic array multi-dimensional virtual domainand cross-correlation tensor statistics is established, a virtual domainsignal tensor is deduced by transformation of cross-correlation tensorstatistics, and cross-correlation virtual domain tensor representationwith a multi-dimensional space information structure retained isconstructed; identifiability of a signal source is ensured byconstructing a signal source identification mechanism for a virtualdomain signal tensor without introducing a spatial smoothing step; andthe Nyquist matched direction-of-arrival estimation is realized by meansof virtual domain signal tensor decomposition.

The above merely shows preferred embodiments of the present disclosure.Although the present disclosure has been disclosed as above in thepreferred embodiments, it is not used to limit thereto. Withoutdeparting from the scope of the technical solution of the presentdisclosure, any person skilled in the art can make many possiblevariations and modifications to the technical solution of the presentdisclosure by using the methods and technical contents disclosed above,or modify it into equivalent embodiments with equivalent changes.Therefore, any simple alteration, equivalent variation and modificationof the above embodiments according to the technical essence of thepresent disclosure without departing from the technical solution of theinvention still fall in the protection scope of the technical solutionof the present disclosure.

1. A three-dimensional co-prime cubic array direction-of-arrivalestimation method based on a cross-correlation tensor, comprising thefollowing steps: (1) constructing at a receiving end with

+

−1 physical antenna array elements in accordance with a structure of athree-dimensional co-prime cubic array, wherein

and

,

and

, and

and

are respectively a pair of co-prime integers, and the three-dimensionalco-prime cubic array is decomposable into two sparse and uniform cubicsubarrays

and

; (2) supposing there are K far-field narrowband non-coherent signalsources from a direction of {(θ₁,φ₁), . . . , (θ_(K),φ_(K))}, carryingout modeling on a receiving signal of the sparse and uniform cubicsubarray

of the three-dimensional co-prime cubic array via a four-dimensionaltensor

∈

×

×

×T (T is a number of sampling snapshots) as follows:

=Σ_(k=1) ^(K)

(μ_(k))∘

(ν_(k))∘

(ω_(k))∘s _(k)+

, wherein, s_(k)=[s_(k,1), s_(k,2), . . . , s_(k,T]) ^(T) is amulti-snapshot sampling signal waveform corresponding to a k^(th)incident signal source, (⋅)^(T) represents a transposition operation, ∘represents an external product of vectors,

∈

×

×

×T is a noise tensor mutually independent from each signal source,

(μ_(k)),

(ν_(k)) and

(ω_(k)) are steering vectors of the three-dimensional sparse and uniformcubic subarray

in an x axis, a y axis and a z axis respectively, and a signal sourcecorresponding to a direction-of-arrival of (θ_(k),φ_(k)) is representedas:a_(x)^((ℚ_(i)))(μ_(k)) = [1, e^(−jπx_(ℚ_(i))⁽²⁾μ_(k)), …, e^(−jπx_(ℚ_(i))^((M_(x)^((ℚ_(i)))))μ_(k))]^(T),a_(y)^((ℚ_(i)))(v_(k)) = [1, e^(−jπy_(ℚ_(i))⁽²⁾v_(k)), …, e^(−jπy_(ℚ_(i))^((M_(y)^((ℚ_(i)))))v_(k))]^(T),a_(z)^((ℚ_(i)))(ω_(k)) = [1, e^(−jπz_(ℚ_(i))⁽²⁾ω_(k)), …, e^(−jπz_(ℚ_(i))^((M_(z)^((ℚ_(i)))))ω_(k))]^(T),wherein,

(i₁=1,2, . . . ,

)

(i₂=1,2, . . . ,

) and

(i₃=1, 2, . . . ,

) respectively represent actual locations of

in i₁ ^(th), i₂ ^(th), and i₃ ^(th) physical antenna array elements inthe x axis, the y axis and the z axis, and

=

=

=0, μ_(k)=sin(φ_(k))cos(θ_(k)), ν_(k)=sin(φ_(k))sin(θ_(k)),ω_(k)=cos(φ_(k)), j=√{square root over (−1)}; (3) based onfour-dimensional receiving signal tensors

and

of the two three-dimensional sparse and uniform cubic subarrays

and

, solving their cross-correlation statistics to obtain a six-dimensionalspace information-covered second-order cross-correlation tensorℛ_(ℚ₁ℚ₂) ∈ ℂ^(M_(x)^((ℚ₁)) × M_(y)^((ℚ₁)) × M_(z)^((ℚ₁)) × M_(x)^((ℚ₂)) × M_(y)^((ℚ₂)) × M_(z)^((ℚ₂))):${\mathcal{R}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}} = {{E\lbrack \langle {\mathcal{X}_{{\mathbb{Q}}_{1}},\mathcal{X}_{{\mathbb{Q}}_{2}}^{*}} \rangle_{4} \rbrack} = {{\sum\limits_{k = 1}^{K}{\sigma_{k}^{2}{{a_{x}^{({\mathbb{Q}}_{1})}( \mu_{k} )} \circ {a_{y}^{({\mathbb{Q}}_{1})}( v_{k} )} \circ {a_{z}^{({\mathbb{Q}}_{1})}( \omega_{k} )} \circ {a_{x}^{{({\mathbb{Q}}_{2})}^{*}}( \mu_{k} )} \circ {a_{y}^{{({\mathbb{Q}}_{2})}^{*}}( v_{k} )} \circ {a_{z}^{{({\mathbb{Q}}_{2})}^{*}}( \omega_{k} )}}}} + \mathcal{N}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}}}}},$wherein, σ_(k) ²=E[s_(k)s_(k)*] represents a power of a k^(th) incidentsignal source,

=E[<

,

>₄] represents a six-dimensional cross-correlation noise tensor,<⋅,⋅>_(r) represents a tensor contraction operation of two tensors alonga r^(th) dimension, E[⋅] represents an operation of taking a mathematicexpectation, and (⋅)* represents a conjugate operation; asix-dimensional tensor

merely has an element with a value of σ_(n) ² in a (1, 1, 1, 1, 1,1)^(th) location, σ_(n) ² representing a noise power, and with a valueof 0 in other locations; (4) as a first dimension and a fourth dimensionof the cross-correlation tensor

represent space information in a direction of the x axis, a seconddimension and a fifth dimension represent space information in adirection of the y axis, and a third dimension and a sixth dimensionrepresent space information in a direction of the z axis, definingdimension sets

₁={1,4},

₂={2,5} and

₃={3,6}, and carrying out tensor transformation of dimension merging onthe cross-correlation tensor

to obtain a virtual domain second-order equivalent signal tensor

∈

×

×

:${\mathcal{U}_{\mathbb{W}}\overset{\bigtriangleup}{=}{\mathcal{R}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2_{\{{{\mathbb{J}}_{1},{\mathbb{J}}_{2},{\mathbb{J}}_{3}}\}}}} = {\sum_{k = 1}^{K}{\sigma_{k}^{2}{{b_{x}( \mu_{k} )} \circ {b_{y}( v_{k} )} \circ {b_{z}( \omega_{k} )}}}}}},$wherein, b_(x)(μ_(k))=

(μ_(k))⊗

(μ_(k)), b_(y)(ν_(k))=

(ν_(k))⊗

(ν_(k)) and b_(z)(ω_(k))=

(ω_(k))⊗

(ω_(k)) respectively construct augmented virtual arrays in thedirections of the x axis, the y axis and the z axis through formingarrays of difference sets on exponential terms, b_(x)(μ_(k)) ,b_(y)(ν_(k)) and b_(z)(ω_(k)) are respectively equivalent to steeringvectors of the virtual arrays in the x axis, the y axis and the z axisto correspond to signal sources in a direction-of-arrival of(θ_(k),φ_(k)), and ⊗ represents a product of Kronecker, so that

corresponds to an augmented three-dimensional virtual non-uniform cubicarray

;

comprises a three-dimensional uniform cubic array

with (3

−

+1)×(3

−

−1)×

−

+1) virtual array elements, represented as:

={(x,y,z)|x=p_(x)d,y=p_(y)d,z=p_(z)d,−

≤p_(x)≤−

+2

,−

≤p_(y)≤−

+2

,−

≤p_(z)≤−

+2

}, The equivalent signal tensor

∈

−

+1)×(3

−

+1)×(3

−

+1) of the three-dimensional uniform cubic array

is modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘b _(y)(ν_(k))∘b _(z)(ω_(k)), wherein,${{{\overset{\_}{b}}_{x}( \mu_{k} )} = {\lbrack {e^{- j{\pi({- M_{x}^{({\mathbb{Q}}_{1})}})}\mu_{k}},e^{- j{\pi({{- M_{x}^{({\mathbb{Q}}_{1})}} + 1})}\mu_{k}},\ldots,e^{- j{\pi({{- M_{x}^{({\mathbb{Q}}_{2})}} + {2M_{x}^{({\mathbb{Q}}_{1})}}})}\mu_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{x}^{({\mathbb{Q}}_{1})} - M_{x}^{({\mathbb{Q}}_{2})}} + 1}}},{{{\overset{\_}{b}}_{y}( v_{k} )} = {{\lbrack {e^{- j{\pi({- M_{y}^{({\mathbb{Q}}_{1})}})}v_{k}},e^{- j{\pi({{- M_{y}^{({\mathbb{Q}}_{1})}} + 1})}v_{k}},\ldots,e^{- j{\pi({{- M_{y}^{({\mathbb{Q}}_{2})}} + {2M_{y}^{({\mathbb{Q}}_{1})}}})}v_{k}}} \rbrack^{T} \in {{\mathbb{C}}^{{3M_{y}^{({\mathbb{Q}}_{1})}} - M_{y}^{({\mathbb{Q}}_{2})} + 1}{and}{{\overset{\_}{b}}_{z}( \omega_{k} )}}} = {\lbrack {e^{- j{\pi({- M_{z}^{({\mathbb{Q}}_{1})}})}\omega_{k}},e^{- j{\pi({{- M_{z}^{({\mathbb{Q}}_{1})}} + 1})}\omega_{k}},\ldots,e^{- j{\pi({{- M_{z}^{({\mathbb{Q}}_{2})}} + {2M_{z}^{({\mathbb{Q}}_{1})}}})}\omega_{k}}} \rbrack^{T} \in {\mathbb{C}}^{{3M_{z}^{({\mathbb{Q}}_{1})}} - M_{z}^{({\mathbb{Q}}_{2})} + 1}}}}$represent steering vectors of the three-dimensional virtual uniformcubic array

in the x axis, the y axis and the z axis corresponding to signal sourcesin the direction-of-arrival of (θ_(k), φ_(k)); (5) as a mirror imageportion

_(sym), of the three-dimensional virtual uniform cubic array

is represented as:

_(sym)={(x,y,z)|x={hacek over (p)} _(x) d,y={hacek over (p)} _(y)d,z={hacek over (p)} _(z) d,

−2

≤{hacek over (p)} _(x)≤

−2

≤{hacek over (p)} _(y)≤

,

−2

≤{hacek over (p)} _(z)≤

}, carrying out transformation by using the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

to obtain an equivalent signal tensor

_(sym)∈

−

+1)×(3

−

+1)×(3

−

+1) of a three-dimensional mirror image virtual uniform cubic array

_(sym), specifically comprising: carrying out a conjugate operation onthe three-dimensional virtual domain signal tensor

to obtain

, carrying out position reversal on elements in the

along directions of three dimensions successively so as to obtain theequivalent signal tensor

_(sym) corresponding to the

_(sym) ; superposing the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

and the equivalent signal tensor

_(sym) of the mirror image virtual uniform cubic array

_(sym) in the fourth dimension to obtain a four-dimensional virtualdomain signal tensor

∈

−

+1)×(3

−

+1)×(3

−

+1)×2, modeled as:

=Σ_(k=1) ^(K)σ_(k) ² b _(x)(μ_(k))∘b _(y)(ν_(k))∘b_(z)(ω_(k))∘c(μ_(k),ν_(k),ω_(k)), wherein,c(μ_(k), v_(k), ω_(k)) = [1, e^(−jπ((M_(x)^((ℚ₂)) − M_(x)^((ℚ₁)))μ_(k) + (M_(y)^((ℚ₂)) − M_(y)^((ℚ₁)))v_(k) + (M_(z)^((ℚ₂)) − M_(z)^((ℚ₁)))ω_(k)))]^(T)is a three-dimensional space mirror image transformation factor vector;(6) carrying out CANDECOMP/PARACFAC decomposition on thefour-dimensional virtual domain signal tensor

to obtain factor vectors b _(x)(μ_(k)),b _(y)(ν_(k)), b _(z)(ω_(k)) andc(μ_(k),ν_(k),ω_(k)), k=1,2, . . . , K, corresponding tofour-dimensional space information, and constructing a signal subspaceV_(s)∈

^(V×K) through a form of their Kronecker products: V_(s)=orth([b_(x)(μ₁)⊗b _(y)(ν₁)⊗b _(z)(ω₁)⊗c(μ₁,ν₁,ω₁), b _(x)(μ₂)⊗b_(y)(ν₂)⊗b(ω₂)⊗c(μ₂,ν₂,ω₂), . . . , b _(x)(μ_(K))⊗b _(y)(ν_(K))⊗b_(z)(ω_(K))⊗c(μ_(K),ν_(K),ω_(K))]), wherein, orth(⋅) represents a matrixorthogonalization operation, V=2(3

−

+1)(3

−

+1)(3

−

+1); by using V_(n)∈

^(V×(V−K)) to represent a noise subspace, V_(n)V_(n) ^(H) is obtained byV_(s):V _(n) V _(n) ^(H) =I−V _(s) V _(s) ^(H), wherein, I represents a unitmatrix; (⋅)^(H) represents a conjugate transposition operation; and (7)traversing a two-dimensional direction-of-arrival of ({tilde over (θ)},{tilde over (φ)}), calculating corresponding parameters {tilde over(μ)}_(k)=sin({tilde over (θ)}_(k)), {tilde over (ν)}_(k)=sin({tilde over(φ)}_(k))sin({tilde over (θ)}_(k)) and {tilde over (ω)}_(k)=cos({tildeover (φ)}_(k)), and constructing a steering vector

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))∈

^(V) corresponding to the three-dimensional virtual uniform cubic array

, represented as:

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))=b_(x)({tilde over (μ)}_(k))⊗b _(y)({tilde over (ν)}_(k)) ⊗b _(z)({tildeover (ω)}_(k))⊗c({tilde over (μ)}_(k), {tilde over (ν)}_(k), {tilde over(ω)}_(k)), wherein, {tilde over (θ)}∈[−90°,π°], {tilde over(φ)}∈[0°,180°], A three-dimensional spatial spectrum

({tilde over (θ)},{tilde over (φ)}) is calculated as follows:

({tilde over (θ)},{tilde over (φ)})=1/(

^(H)({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))(V_(n) V _(n) ^(H))

({tilde over (μ)}_(k),{tilde over (ν)}_(k),{tilde over (ω)}_(k))),Spectral peak search is carried out on the three-dimensional spatialspectrum

({tilde over (θ)},{tilde over (φ)}) to obtain a direction-of-arrivalestimation result.
 2. The three-dimensional co-prime cubic arraydirection-of-arrival estimation method based on a cross-correlationtensor according to claim 1, wherein the structure of thethree-dimensional co-prime cubic array in step (1) is described as: apair of three-dimensional sparse and uniform cubic subarrays

and

are constructed in a rectangular coordinate system, wherein

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, they axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(1x),

dm_(1y),

dm_(1z)), m_(1x)=0,1, . . . ,

−1, m_(1y)=0,1, . . . ,

−1, m_(1z)−0,1, . . . ,

−1};

comprises

×

×

antenna array elements, with array element spacings in the directions ofthe x axis, the y axis and the z axis being

d,

d and

d respectively, with locations in the rectangular coordinate systembeing {(

dm_(2x),

dm_(2y),

dm_(2z)), m_(2x)=0,1, . . . ,

−1, m_(2y)=0,1, . . . ,

−1, m_(2z)=0,1, . . . ,

−1}; a unit spacing d has a value half of an incident narrowband signalwavelength) λ, i.e., d=λ/2; subarray combination is carried out on the

and

in such a way that array elements on the (0, 0, 0) location in therectangular coordinate system are overlapped so as to obtain athree-dimensional co-prime cubic array actually containing

+

−1 physical antenna array elements.
 3. The three-dimensional co-primecubic array direction-of-arrival estimation method based on across-correlation tensor according to claim 1, wherein the second-ordercross-correlation tensor statistics of the three-dimensional co-primecubic array in step (3) are estimated by calculating cross-correlationstatistics of T sampling snapshots of the receiving signal tensors

(t) and

(t) in reality:${\hat{\mathcal{R}}}_{{\mathbb{Q}}_{1}{\mathbb{Q}}_{2}} = {\frac{1}{T}{\sum_{t = 1}^{T}{{{\mathcal{X}_{{\mathbb{Q}}_{1}}(t)} \circ {\mathcal{X}_{{\mathbb{Q}}_{2}}^{*}(t)}}.}}}$4. The three-dimensional co-prime cubic array direction-of-arrivalestimation method based on a cross-correlation tensor according to claim1, wherein the equivalent signal tensor

of the three-dimensional virtual uniform cubic array

in step (4) can be obtained by selecting elements in the equivalentsignal tensor

of the three-dimensional virtual nonuniform cubic array

corresponding to locations of virtual array elements in the

.
 5. The three-dimensional co-prime cubic array direction-of-arrivalestimation method based on a cross-correlation tensor according to claim1, wherein in step (6), CANDECOMP/PARAFAC decomposition is carried outon the four-dimensional virtual domain signal tensor

to obtain factor matrixes B _(x)=[b _(x)(μ₁), b _(x)(μ₂), . . . b_(x)(μ_(K))], B _(y)=[b _(y)(ν₁), b _(y)(ν_(z)), . . . , b _(y)(ν_(K))],B _(z)=[b _(z)(ω₁), b _(z)(ω₂), . . . , b _(z)(ω_(K))] andC=[c(μ₁,ν₁,ω₁), c(μ₂,ν₂,ω₂), . . . , c(μ_(K),ν_(K),ω_(K))], wherein,CANDECOMP/PARAFAC decomposition of the four-dimensional virtual domainsignal tensor

follows a uniqueness condition as follows:

_(rank)(B _(x))+

_(rank)(B _(y))+

_(rank)(B _(z))+

_(rank)(C)≥2K+3, wherein,

_(rank)(⋅) represents a Kruskal rank of a matrix, and

_(rank)(B _(x))=min (3

−

+1, K),

_(rank)(B _(y))=min(3

−M

−1, K),

_(rank)(B _(z))=min(3

−

+1, K),

_(rank)(C)=min(2, K), min (⋅) represents an operation of taking aminimum value; when spatial smoothing is not introduced to process thededuced four-dimensional virtual domain signal tensor

, a uniqueness inequation of the above CANDECOMP/PARACFAC decompositionis established.
 6. The three-dimensional co-prime cubic arraydirection-of-arrival estimation method based on a cross-correlationtensor according to claim 1, wherein in step (7), a process of obtaininga direction-of-arrival estimation result by three-dimensional spatialspectrum search specifically comprises: fixing a value of {tilde over(φ)} at 0°, gradually increasing {tilde over (θ)} to 90° from −90° at aninterval of 0.1°, increasing the {tilde over (θ)} to 0.1° from 0°,increasing the θ to 90° from −90° at an interval of 0.1° once again, andrepeating this process until the {tilde over (φ)} is increased to 180°,calculating a corresponding

({tilde over (θ)},{tilde over (φ)}) in each two-dimensionaldirection-of-arrival of ({tilde over (θ)},{tilde over (φ)}) so as toconstruct a three-dimensional spatial spectrum on a two-dimensionaldirection-of-arrival plane; searching peak values of thethree-dimensional spatial spectrum

({tilde over (θ)},{tilde over (φ)}) in the two-dimensionaldirection-of-arrival plane, permutating response values corresponding tothese peak values in a descending order, and taking two-dimensionalangle values corresponding to first K spectral peaks as thedirection-of-arrival estimation result of a corresponding signal source.